library(lmtest)
library(sandwich)
library(foreign)



# Load data ---------------------------------------------------------------

svyww <- read.csv("../data/asylum_data.csv")


# Country vector ----------------------------------------------------------

countries<-c("Germany","Hungary","Sweden","Austria","Norway",
             "Switzerland","Denmark","Netherlands","Greece","Czech Republic","Italy",
             "Poland","Spain","France","United Kingdom")

countries <- sort(countries)

vars <- c("Female","Age","Employed","EISCED","IncomeDecile","IdeoScale","Knowledge")

svyww$tcondition <- NA
svyww$tcondition[svyww$InfoTrt == 0 & svyww$ConsTrt == 0] <- "A"
svyww$tcondition[svyww$InfoTrt == 1 & svyww$ConsTrt == 0] <- "B"
svyww$tcondition[svyww$InfoTrt == 0 & svyww$ConsTrt == 1] <- "C"
svyww$tcondition[svyww$InfoTrt == 1 & svyww$ConsTrt == 1] <- "D"

bal.list <- list()
for (i in 1:length(countries)){
  
  svycty <- subset(svyww,svyww$cty == countries[i])
  df <- matrix(NA, nrow = length(vars)*2 + 2, ncol = 5, dimnames = list(c("N","empty",rep(vars,each=2)),c("A","B","C","D","ANOVA p val")))
  
  for (j in 1:length(vars)){
    
    df[1,"A"] <- sum(svycty$tcondition == "A")
    df[1,"B"] <- sum(svycty$tcondition == "B")
    df[1,"C"] <- sum(svycty$tcondition == "C")
    df[1,"D"] <- sum(svycty$tcondition == "D")
    
    k <- j*2 + 1
    
    df[k,"A"] <- mean(svycty[svycty$tcondition == "A",vars[j]],na.rm=T)
    df[k+1,"A"] <- sd(svycty[svycty$tcondition == "A",vars[j]],na.rm=T)/
      sqrt(sum(!is.na(svycty[svycty$tcondition == "A",vars[j]])))
    
    df[k,"B"] <- mean(svycty[svycty$tcondition == "B",vars[j]],na.rm=T)
    df[k+1,"B"] <- sd(svycty[svycty$tcondition == "B",vars[j]],na.rm=T)/
      sqrt(sum(!is.na(svycty[svycty$tcondition == "B",vars[j]])))
    
    df[k,"C"] <- mean(svycty[svycty$tcondition == "C",vars[j]],na.rm=T)
    df[k+1,"C"] <- sd(svycty[svycty$tcondition == "C",vars[j]],na.rm=T)/
      sqrt(sum(!is.na(svycty[svycty$tcondition == "C",vars[j]])))
    
    df[k,"D"] <- mean(svycty[svycty$tcondition == "D",vars[j]],na.rm=T)
    df[k+1,"D"] <- sd(svycty[svycty$tcondition == "D",vars[j]],na.rm=T)/
      sqrt(sum(!is.na(svycty[svycty$tcondition == "D",vars[j]])))
    
    mod <- lm(svycty[,vars[j]] ~ svycty[,"tcondition"])
    df[k,"ANOVA p val"] <- anova(mod)$`Pr(>F)`[1]
    
  }
  
  df <- round(df,2)

  rownames(df) <- NULL
  
  df <- cbind(Variable = c("N","empty",rep(vars,each=2)),as.data.frame(df))
  df$Variable <- as.character(df$Variable)
  
  for (j in 1:length(vars)){
    df$Variable[2*j+2] <- " "
    df[2*j+2,2:5]  <- paste("(",df[2*j+2,2:5],")",sep="")
  }
  df$Variable[2] <- " "
  
  names(df) <- c("Variable","Condition A","Condition B","Condition C","Condition D","ANOVA pval")
  df$Variable <- c("N","","Female","","Age","","Employed","","Education","","Income Decile","","Political Ideology","","Knowledge Index","")
  
  bal.list[[i]] <- df
  names(bal.list)[[i]] <- countries[i]
  
}


for (j in countries){
  print(j)
  print(bal.list[[j]])
  cat("\n\n")
}

